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1. INTRODUCTION 


Among the various types of available data relevant to the establishment 
of geometric control on the moon, the only one covering significant portions 
of the lunar surface (20 %) with sufficient information content, is lunar photog- 
raphy, taken at the proximity of the moon from lunar orbiters . The naturally 
suited tool for the analysis of such data is the well known method of analytical 
phototriangulation. However, analytical phototriangulation in its earth-bound 
applications has not been traditionally viewed by photogrammetrists as a method 
to establish, but rather to extend or density already existing geodetic control. 
The lack of high quality geometric control on the moon-requires the introduc- 
tion of orbital control which can be established from earth-based tracking of 
lunar orbiters using models containing information on the lunar ephemeris, 
on the rotations of the earth and the moon, station positions on the earth and 
on the moon and on the lunar gravity field. 

Past experience from lunar phototriangulation with separately deter- 
mined and constrained orbital control shows unacceptable inconsistencies [7], 
and there has been a call for the employment of advanced orbit determination 
techniques emphasizing the geometric aspects rather than the previously em- 
ployed gravimetric aspects. Such techniques have become available mainly 
through the philosophy expressed in the work of Ingram and Tapley [9] . 

This work attempts to explore certain new ideas aimed in somewhat 
different directions . The systematic pattern in some of the inconsistencies 
mentioned above [7, Fig. 5 and 6] • strongly indicates that they are at least 
partly due to referencing the orbit geometry to poorly determined frames, 
e. g. , to the principal moments of inertia axes or to the selenographic systems . 

In this paper the idea of "free geodetic networks" [4] is introduced as 
a tool for the statistical comparison of the geometric aspects of the various 
data used. Further, methods are developed for the updating of the statistics 
of observations and the a priori parameter estimates to obtain statistically 
consistent solutions by means of the "optimum relative weighting" concept. 
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2 . REVIE W O F A VAIL AB LE DATA 

/ 

/ 

Although other types of data such as laser ranging to the moon and differ- 
ential VLBI are relevant to the establishment of selenodetic control, the empha- 

1 5 

sis here is on the utilization of orbital metric photography in conjunction with 
supporting laser altimetry and stellar photography and of a priori estimates of 
satellite positions. More specifically, the data here considered have been ob- 
tained by the Mapping Camera System aboard the spacecraft during the Apollo 
15, 16 and 17 missions. Since a detailed description of the system can be found 
in [13] and [3], only a general description is given here. 

The Mapping Camera System is composed of a terrain camera, an associ- 
ated stellar camera and a laser altimeter. The terrain camera has a focal length 
of 76 mm, a 115 mm x 115 mm format and reseau marks with a spacing of 10 mm 
recorded on every photograph for film deformation compensation. The camera 
is capable of compensating for forward motion and the shutter interval is auto- 
matically set by means of separate sensor measuring the brightness of the lunar 
surface. Resolution depends on the solar altitude with an optimum of 90 lines 
per mm. The stellar camera has a 76 mm focal length, a format of 32mm X 
25mm and a reseau grid with 5mm spacing. The laser altimeter records the 
spacecraft altitude with an accuracy of ± 2m. 

During each revolution as the spacecraft passes over the sun -illuminated 
part of the lunar surface, the terrain and stellar cameras and the altimeter oper- 
ate simultaneously to provide three types of observations at each exposure time. 
The terrain camera takes a strip of almost vertical lunar photographs , while the 
stellar camera takes a corresponding series of star photographs in the direction 
of the flight about 4° to 6° above the horizon. At the same time the altimeter 
measures the distance between the camera exposure point and the intersection 
of the camera axis with the lunar surface. The terrain camera also records 
this laser illuminated point. From the star field photographs, the orientation 
of the stellar camera is determined with respect to a "star catalog" system and 
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the orientation of the terrain camera in the same system can be determined 
using the relative orientation of the two cameras known from preflight cali- 
brations . 

Earth-based tracking of the spacecraft is used to determine the orbit 
geometry, more specifically, the coordinates of the camera with respect to 
some moon-fixed coordinate system at each exposure time. 
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3. COMMENTS ON COORDINATE SYSTEMS 
3. 1 Coordinate Systems in General 

The concept of geodetic control traditionally has always been connected 
to the use of coordinates with respect to some frame of reference. Although 
coordinates are not the only necessary means to represent geodetic control, 
they have been introduced in practice as a matter of convenience because they 
allow the use of Cartesian analytic geometry, inspite of possible pitfalls. 

Reference frames can be divided into natural and conventional frames . The 
former naturally arise from the physics of a given situation (e.g. , principal 
axes of inertia), while the latter are arbitrarily set to meet some criteria. 
Coordinates have a physical meaning and they can be determined (estimated) 
from observations only when they are referred to the natural frames of 
reference. 

While some observations are invariant under coordinate system trans- 
formations and therefore they can be analyzed in any arbitrary reference frame, 
other observations are connected with certain reference frames and must be 
analyzed therein. For example, while range or range rate observations are 
coordinate system invariant, the stellar camera observations are naturally 
variant with respect to the "star catalog" frame of reference. 

Modeling the physical processes is also connected to certain reference 
frames. For example, lunar theory can be viewed as the time history of 
the geocenter to the lunar center of mass vector, and the physical librations 
of the moon are connected to the principal axes of inertia selenocentric ref- 
erence frame. For the foregoing reasons it seems advisable to define precisely 
the reference frames to be used in this report. 

3. 2 Traditional Coordinate Systems 

The first reference frames for the description of coordinates of lunar features 
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emerged at times when the only possible means of observations were the low 
accuracy earth-based optical observations, such as heliometry and later lunar 
photography. In view of the low accuracy of these observations, some aspects 
in the precise definition of coordinate systems may have seemed to be too theore- 
tical and/or irrelevant. As a result, errors are frequent in selenographic papers 
and the confusion found its way even into astronomical almanacs [6]. The main 
problem arises from the failure to distinguish between the following four direc- 
tions (supposedly) through the selenocenter: 

1) The direction of the axis of maximum moment of inertia. 

2) The direction of the instantaneous axis of lunar rotation. 

3) The Cassini axis, i.e., the axis whose change of orientation with 
respect to an inertial (ecliptic) system is governed by the three 
empirical laws of Cassini. 

4) The direction of the mean lunar rotation axis . 

A detailed description of the current and past situation is given in a paper by 
Habibullen (1971) [6]. His terminology is adopted also for this paper. The 
following four selenocentric systems are considered: 

1) True Selenographic System — The z axis is the instantaneous rotation 
axis of the moon. 

2) Dynamical System — This one is usually referred to as the principal 
axes of inertia system. It is defined as the Cartesian selenocentric system - 
for which the following products of inertia vanish: 

D = / M yz dm = 0, E = Xi xz dm =0, F = fx y dm = 0. 

3) Cassini Selenographic System — This system is such that its rotation 
with respect to the ecliptic system is given at any epoch by three Eulerian angles 
defined as follows: 

if)o = Q 

(Do * 180°+ (£ J + Q) 

8o = I 

where Q is the mean longitude of the ascending node of the lunar orbit; J the 
mean lunar longitude and I the inclination of the Cassini equator (not the orbital 
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plane) with respect to the ecliptic. 

4) Me«n S el eno graphic System — The z axis is the mean rotation of the moon. 

For systems (2) and (3) the direction of the x axis is defined. For systems (1) and 
(4) the "zero point" on the corresponding equator remains to be defined or to be arbi- 
trarily selected. The only moon-fixed coordinate system by definition is the 
Dynamical System. Systems (1) and (3) can also be considered moon-fixed if 
they are referred to a certain reference epoch. System (4) has the most compli- 
cated definition and though it is the most widely used, it is the least well realized 
(established). Formally, a mean position is to be defined by means of a time 
average of the time variant positions, i. e.,, it involves an integral with respect 
to time over the interval ( - ro , +“•).. Such an integration could be carried out 
only if the variation of the position of the instantaneous axis of rotation with respect 
to the rigid lunar body has a purely periodic character. Since this is not neces- 
sarily the case, the integration has to be realized over a certain time interval 
and then system (4) becomes moon -fixed with respect to this certain fixed time 
interval . 

The Dynamical (coordinate) System is related to the Cassini Selenographic 
System through the relations [ 6 ]: 

= 0o + O' 

<P = <Po + (r - cr) 

9 = d 0 + p 

where cr, p and r are the physical librations in node, inclination and longitude, 
respectively. The difficulty in detecting the position of the principal axes of 
inertia by means of present observations -necessitates the use of an "observable" 
system such as the Mean Selenographic or the True Selenographic Systems 
referred to some fixed epoch. In practice, therefore, the physical libration 
parameters cr , p and- T connect the Cassini System to one of the above systems. 

In earliest works the mean rotation- axis of the moon (Mean Selenographic 
System) was taken as identical to the principal axis of maximum inertia (Dynamic 
System) [6], The distinction between these systems is not "critical when construc- 
ting a libration theory in which the product of inertia D, E, F are to be constrained 
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to zero, because the products of inertia in the Mean Selenographic System are 
close enough to zero. The detection of the directions of the principal axes of 
inertia with respect to another moon-fixed system is possible only through the 
knowledge of the gravity field, specifically the second degree harmonics, with 
respect to the moon-fixed system. In Appendix C the possibility and methods 
of determining these directions are dealt with separately. 

3. 3 Natural and Arbitrary Coordinate Systems 

From the point of view of establishing geometric control, the primary 
objective is to determine the shape and scale of a network formed by a cluste: 
of points on the lunar surface. The shape and scale of a network is determined 
when the following quantities can be determined: 

1) The angLe Pj Pj P k between lines P t Pj , Pj P k through any points 
Pi , Pj , P k in the network. 

2) The distance d t j of the line segment PjPj j oining any two points 
Pi and Pj in the network. 

Group (X) quantities determine only the shape of the network. 

Alternatively, for representing the network geometry one may use the 
coordinates of the points with respect to a coordinate system which is invariant 
in space with respect to the rigid network of points. In this case both the angles 

A 

the distances d t j are functions of the Cartesian coordinates of the 
points, which functions are invariant with respect to the choice of a particular 
coordinate system. This reflects the fact that the shape and scale of the network 
is independent of the choice of the reference system. 

The coordinates of the points with respect to some particular reference 
frame contain more information than specified under shape and scale. They 
also define the position of the coordinate system with respect to the network of 
points. If the coordinate system is a natural one, i.e. , if it is the consequence 
of the natural characteristics of the physical objects involved in the model, and 
the problem is to determine the position of a network of points with respect to 
this system, one should proceed in two steps, as follows: 
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1) Find the shape and scale of the network of points . 

2) Find the position of the natural system of reference with respect 
to the network! 

A minimum number of angles and distances which uniquely determine the shape 
and scale of the network is called a "fundamental set." All other angles and 
distances can be determined as functions of the fundamental set. Such a funda- 
mental set also serves as a representation of the shape and scale. If coordi- 
nates are to be used as a representation of the shape and scale only, the refer- 
ence frame’s position with respect to the network of points must be known by 

definition. Such a reference frame can be established by means of the concept 

* 

of minimal constaints on the coordinates . A set of coordinate conditions which 
specify the positions of the reference frame with respect to the network is 
called "minimal constraints. " An' easy way to set such a constraint is to 
assign constant values to six coordinates distributed ober at least three differ- 
ent points. For a more detailed discussion on minimal constraints see [4], 

[2] and [14]. 

3.4 Role of Coordinates in Selenodesy 

The problem of finding the position of a network of points with respect 
to a natural system can now be viewed in two parts: 

1) Find the coordinates of the points in the network with respect to 

an arbitrary reference system introduced through the use of a set of minimal 
constraints . 

2) Find a set of transformation parameters which specifies the relative 
positions of the natural coordinate system with respect to the arbitrary system. 
The natural choice of such parameters is of course, three translations and 
three orientation angles. 

The reason for the above dichotomy becomes obvious if one brings to 
mind that beyond the deterministic definitions concerning positions of points, 
one in practice has to extract estimates of these positions from observations 
which can be either coordinate system dependent or not. Furthermore, 
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observations may depend on more than one reference system. The situation 
becomes even more complicated when the relative position of these systems 
varies with time. For example, the representation of the geometry of a system 
of points on the moon requires a moon-fixed system of reference. Such a system 
can be either arbitrary (defined by means of a set of minimal constraints), or 
natural. The only natural system qualifying is the Dynamical System. Its 
advantage is that its position depends only on distribution of masses and that 
it is time invariant. Lunar photography and laser altimetry provide us with 
data which are independent of coordinate systems and, therefore, their analysis 
could be performed in an arbitrary system. However, because stellar photography 
is star catalog system dependent and the libration theory transforms it to a moon- 
fixed system, the lunar photography and laser data could also be analyzed in this 
system. Estimates of orbital geometry are given in the form of camera expo- 
sure station positions with respect to some moon-fixed system, not necessarily 
the same used in the libration theory for the transformation of stellar camera 
observations. 

For the purpose of estimating only the shape and scale of the network 
it is quite easy to ’’free" the observational data from their dependence on 
coordinate systems. The problem, however, of relating the network to any 
of the data reference frames still remains. 

3.5 Coordinate Systems in the External Information Used 

At this point it is appropriate to say a word about the theories of lunar 
motion and physical libration. A lunar theory is a solution to the problem of 
the variation of the position (geocenter-selenocenter) vector with time . As 
such, it is an integral in one form or another of the differential equations 
governing the motion of this vector-. The solution contains some integration 
constants which are either the values of the initial state of the Dynamical 
System (numerical integration) or they are- implicit functions of the initial 
state (analytical theories) . A libration theory is similarly a solution to the 
problem of the rotation of the moon (or more precisely, of a moon-fixed 
system, usually the Dynamical System) about the selenocenter. - 
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In both theories in addition to the problem of determining the positions 
vector or the set of libration angles at any epoch, in terms of the integration 
constants, the values of the integration constants themselves also have to be 
determined. This is possible only through an estimation procedure based upon 
actua l observations. The observations are of a geometric nature, distances 
and angles related to points on the earth and on the moon. The position of 
these points is represented through their coordinates with respect to an earth- 
fixed or moon-fixed system, and a theory (model) must be available for the 
description of the relative motion of the two systems . This model can be sepa- 
rated into three components all with respect to an inertial system: rotation of 
a geocentric earth-fixed system, motion of the geocenter-selenocenter vector 
and the rotation of the moon-fixed selenoeentric system. The choice of the 
geocenter and selenocenter as system origins is natural because lunar theory 
is developed upon differential equations governing the motion of those two 
centers of mass . 

The choice of the Dynamical System among all possible selenoeentric moon- 
fixed system is attractive to libration theorists [12] because of the relatively 
simple form that Euler's dynamical differential equations take when this 
system is used. 

In the process of estimating the integration constants of a lunar or 
libration theory, the coordinates of observational or observed points can be 
estimated as a by-product. Their accuracy will reflect how well the position 
of the coordinate system with respect to these points is known. 

The user of lunar ephemeris or libration theory must realize that the 
position of the selenocenter and the orientation of the principal moments of 
inertia axes given in those theories are only estimates of the position of the 
true selenocenter and the true orientation of the axes. The accuracy of 
those estimates is supposedly given by means of the associated statistics 
(variances and covariances). These statistics, however, are realistic only 
when the model is perfect and there are no computational (round-off) errors 
and the observational input statistics are true. 

Another problem arises from the fact that before a solution is attempted. 
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there is a need for two sources of a priori information: 

1) A libration theory for the reduction of stellar camera observations 
from an inertial ("star catalog") system to a moon-fixed system. 

2) Orbital geometry estimates for the strengthening of the photogram- 
metric solution. 

In general such information is available- from external sources and it 
may or may not' fit the model. In this case either the model needs to be 
changed (but not if it is believed to be correct), or the statistics included 
with the information need to be modified in some appropriate way. 
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4. APPROXIMATE ERROR ANALYSIS OF THE CONTRIBUTION 
OF THE VARIOUS DATA TYPES 

The purpose of this section is to give a.projection on how the different 
types of data contribute to the determination of geometric control on the moon. 
The analysis will be with respect to a strep of photographs covering an arc 
of appro xim ately 180° on the lunar surface. The analysis requires four systems 
of reference: 

System (1) 

An arbitrary system introduced by the following set of minimal constraints 

Xi — Yj. — Z x = Xi = — 00 1 = 0 

where X x , Y 1 , Z x are the coordinates of the first camera exposure 
station and Xi , (fi\ , Otty are the camera orientation angles in the tradi- 
tional photogramme trie sense (X axis positive in the direction of the 
flight, Y perpendicular to the orbital plane and Z in the direction of the 
selenocentric vector pointing away from the moon, forming a right 
handed system). 

System (2) 

The moon-fixed system whose orientation with respect to the ecliptic 
system is given by the physical lib-ration parameters to be used in the 
transformation of stellar camera observations (Dynamical System) . 

System (3) 

The moon-fixed system in which estimates of camera exposure station 
positions are given from a previous orbit analysis. 

System (4) 

A special cylindric type system designed to match the geometry of a 
photo strip. In this system after a point is orthogonally projected into 
the orbital plane, itsX coordinate is the angle between the projection's 
selenocentric vector and the Z axis of system (1) . The Z axis is along 
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Assuming that the only source of errors is the error 6p in the measurements 
of photo coordinates from Figure 2, it is clear that the error in direction X is 

6X = 2 “p 5p . 

i 


Similarly, from the YZ plane 

6 Y = 2 5p. 

t 

From Figure 3 the error components from the same displacements 6p in 
the direction Z are 


6Zi 



H 

Bi 


and 


6Z 2 = H 

6P-|- B 2 


Assuming that B x = B 2 



the total error is 


6 Z = 6Zi + 6Z 2 

az , _L 


From Figure 4 the following obvious relationships can be derived: 


Vh 2 + B 2 6<p! = AXsinG 


sin 0 = 


H 


Vh 2 + B 2 


AX = -|-5p. 


Substituting the last two equations into the first one yields 


H' 


6<Pi H 2 + B 2 


1 

y 6P » 


(4.1) 


(4.2) 


(4.3) 
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Figure 2 illustrates in the XZ plane of system (4), the geometry of two 
vertical photographs taken from the same- height H at a distance B apart. 
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therefore, the errors in the orientation angles <p and oo are 


6<p = 26<Pi = — 


f 2 H 2 5p = 6w . 
f H 2 +B 2 


(4.4) 


From Figure 5 the error in the orientation angle x can be derived as 
follows: 


and thus 


B6x x = AY =~6p,- 
2 H 

6 x = 26xi = ~ — 6p . 


(4.5) 


In accordance with the definition in coordinate system (4), the errors in 
the coordinates Y and Z, contributed by 6x and 6 <p will be 

6Y X = B • 6x (4.6) 


SZp = B • 6<p (4.7) 

As a numerical example close to reality, let Figure 6 show a strip 
of photographs of semicircular shape over an arc of 180°. 


Orbit 
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The height is H - 110 km; the lunar radius is R = 1738km and there is a 
total of 70 photographs in the strip so that 83km. The focal length of 
the camera is 76 mm and the error of plate coordinate measurements is 
6p ~ 5um. The problem is to find the total accumulated errors in AX = 
S6X, etc . at the end point of the strip. 

For two successive photographs, equations 4.1 through 4.3 yield 

6X = 5Y - 14.47m 

6Z ' = 38.36m 

and equations 4.4 through 4.7 give 

OOJ = 6<p = 17?29 

6 X = 17?98 

6Z (0 = 6.96 m 

6 Y„ = 7.24 m. 

/t 

Thus for photography alone, the total accumulated errors at the end point 
(5805 km away) are: 


, 70 

AX = ' £ 6X = 

i 

1.013 km 

or 

174 PPM 

70 

AY - y (6 Y + 

j 

6Y X ) = 1.520 km 

or 

262 ppM 

7C 

r — : 

• 



AZ . = > (6 Z + 

dZ^) •= 3.170 km 

or 

546 ppM. 


i 


The effect of altimetry is mainly in the Z direction . Since the altimetry 
error is 6a =■ 2m, the error in Z for one pair of photographs becomes 
only 

6z = 2 6° = 4 m, 
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and, therefore, the total error in Z reduces to 

AZ = £{6z+ 6Z<p) = 767.13m or 132ppM 

which is considerably smaller than the AZ from photography alone . 

The stellar photography with an accuracy of 20" , does not permit the 
accumulation of errors in 6x, 6cp , 6o), thus 6Y X = SZ^ - 0, and therefore, 
the total errors in Y and Z, including the effects of photography and altimetry, 
will be 

AY = E6Y = 1018.15m or 174ppM 
AZ = L6Z = 280.0 mor 48ppM 

again a considerable reduction. The total error in X will remain the same 

AX = 1013.15m or 174ppM. 

A displacement in X or Y, however, causes a displacement in orientation 
which, with respect to the selenocenter in the XZ plane is 


1.013 km 
1848. km 


11 3 "06 


and one of similar magnitude in the YZ plane. Such angular displacements are 
controlled by the stellar camera orientation and therefore the errors in X or Y 
are limited to 20" x (R+ H), i.e., 

AX = AY = 180m or 31ppM. 

These numbers together with the maximum expected AZ = 280m or 48ppM, 
represent reasonable estimates of the accumulated errors over the strip coming 
from a 5fzm photo-measure error in each photograph. 
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5 . REDUCTION AND MODELING OF THE OBSERVATIONS 
AND THE USE OF ORBITAL INFORMATION 


5.1 The Use of Observations from the Mapping Camera System 


Since the terrain camera is separately calibrated, the coordinates used in 
the analysis are not the measured photo -images coordinates, but those which 
have been corrected for distortion, film shrinkage, etc., i.e., they are assumed 
to be free of systematic errors and as such they fit the simple photogrammetric 
model of the following projection equations: 


= Xjj - Xo + f 


fy = Jii ~ YO + f 

S„ 


In the above equations 


Pi/ 


"AXjj" 

Qu 

= Mj 


Sir 




( 5 . 1 ) 


where Mj = R 3 (Xj )R 2 (<Pj) Ri (tOj > 


AXij” 

l! 

x, - 

i 

AYy 


Y, - 

Yj 

AZjj 
U -1 


. Zl “ 

J 

•^3 

1SJ 


X t , Yj , Z i are the coordinates of i th ground point, 

X j} Yj, Z j , co j , <p j > Xj are the coordinates and orientation angles of j th 
camera exposure station. 


20 



x i j > y i i are the plate coordinates and 

x 0 , yo , f are the plate coordinates of the principal point and the focal length, 
all determined from calibration . 

The altimeter observations are modeled in a straightforward way, 
as follows: 

- [ (X j - Xj) 2 + (Y 1 -Yj) 2 + (Z 1 -Zj) 2 ]^~ 0 (5.2) 

where d ^ is the measured distance between the j th camera exposure station 
and the corresponding (illuminated) i th point on the ground. 

The stellar camera produces observations with respect to a star catalog 
system, which subsequently will be reduced to a moon- fixed system as shown 
below. Using a libration theory we have 


X cj 

X, 

Xs 

Kj , ' , C7 j 


> <Pj » Wj 


We have 


vector in the camera axes system at the j th exposure 
same vector in the inertial system 

same vector in selenocentric moon-fixed system used in 
the libration theory 

camera orientation angles with respect to the inertial 
system at the j 1 h exposure 

same angles with respect to the selenocentric moon-fixed 
system . 


Xcj = M(K j , , 0,)X, 

Xcj = M(Xj , (Pi , C0j)X s . 


At epoch t j the three Euler ian angles from the libration theory are 


and 

or, 

and 


= f(tt) > 9j = ®(t j) , = 0(t j) 

x s = iM^j ) R i (- 9 j)R 3 (?j )x i 

X s = R(0j , 0J , €j )X , 

Xcj = M (Xj „.<pj , COj)R(^j , 0j , €j)X| . 
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Also 


M(K } , , Q,) 


and 


M (Xj , <pj , Wj ) R (0 i , 0! 


*j) 


M (Xj , wj = M (Kj , , a J )R 3 (-ej)Ri(0 i )Rs(-0j ) - 

= M (Kj , , Oj)5 (0 j , 0j , €. ). 


After the matrix M has been computed, the angles Xj , (pj , Wj can be calculated 
from the elements of M as follows: 



The variance-covariance matrix of the quantities Xj , <p } , Wj , as computed 
from the variances and covariances of Ky, $ s , , 0j , €j given in 

Appendix B . 

5. 2 General Comments on the Use of Orbital Information 

The analysis of tracking data for orbit determination depends heavily on the 
gravity field model used. The reference system naturally suited for such an 
orbit analysis is a selenocentric inertial system. The analysis of earth-based 
tracking data involves the rotation of the earth, earth station coordinates and the 
lunar theory. The role of libration theory is limited to the description of varia- 
tions in the moon's gravity field due to the moon's rotation with respect to the 
inertial selenocentric system. 

The geometry of the observations (range, range rate) is poor since tracking 
takes place always from the same direction and it is interrupted while the lunar 
satellite is on the far side of the moon. This interruption of tracking and the fact 
that the knowledge of gravity information on the far side is still rather poor. 
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impose strong limitations on the use of long arcs of more than one revolution'. 

For the purpose of this analysis, side overlaps of photo-strips from dif- 
ferent passes provide better quality ties between the network points than one 
could hope to get through long arc techniques. Each orbital arc corresponding 
to one photo-strip is viewed as a separate orbit and tracking provides informa- 
tion on the geometry of camera points along each arc. 

• Both range and range rate observations can provide relatively good estimates 
on the shape and scale of each arc, but not on rotational displacements with the 
earth (geocenter or tracking station) as a center. Positioning with respect to 
the selenocenter and orientation mainly come from the lunar gravity model. 

The weakness of the orbit in positioning and especially in orientation may well 
not be reflected in the statistics (variances and covariances) of orbit point coor- 
dinates. The reason is that statistics strongLy depend on the model, and they 
are realistic only as far as the model is realistic . While the shape and scale 
of the orbit mostly depend on the observations themselves for which realistic 
a priori statistics are available; the truncation of the gravity field and its 
uncertainties make the statistics of gravity field dependent parameters (position 
and orientation) unrealistic. Past phototriangulation solutions with use of orbital 
support indeed.. show displacements in the "adjusted" orbit, mainly in rotations of 
Ihe-arcs -around some point in the orbit (see e.g. , Figures 5 and 6 in [7] ). 

The conclusion is that in the use of orbital, support in phototriangulation, 
one needs to consider two deviations from the classical implementation of a 
priori estimates and their statistics. First, different passes will have to be 
freed from their reference to a common moon-fixed system. Second, since 
statistics of orbit positions are not reliable, they need to be updated in a proper 
way. In simple terms , the orbit information has to be properly weighted rela- 
tive to the photographic information. In the next chapters these two problems 
are investigated. 
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6. STATISTICAL COMPARISON OF DATA 

6. 1 The Input Data, Their Contributions and Their Associated Reference Systems 

Table 1 is a summary of the various data types and where they contribute to 
the problem of determination of shape, scale, orientation and positioning of a 
selenodetie network. 


Table 1 


Type of Data 

Shape 

Scale 

Orientation 

Positioning 

Lunar Photography 

YES 




Stellar Photography 

YES 


YES 


Altimetry 

YES 

YES 



Orbit State Vectors 

YES 

YES 

YES 

YES 


Since all data types contribute to the shape of the network, it becomes a possible 
tool for detecting inconsistencies between the various types of data. Positioning 
depends only on orbital support and one can expect the tie between the selenodetie 
network and the selenocenter to be only as strong as provided by the orbit analysis. 
Orientation depends both on stellar photography and on the orbit, but the low quality 
of stellar photography makes the detection of possible inconsistencies in orbit 
orientation with respect to a moon-fixed system almost impossible. 

In this work the target of the analysis is only the determination of'the shape 
and scale of the network, and for positioning and orientation, one has to rely on 
orbital support. The latter also provides a means for scale information, and as 
such it can detect scale inconsistencies in the altimetry data. 

To avoid inconsistencies in the orientation and positioning between the different 
passes and also between these passes and the corresponding stellar photography, 
in solutions directed to determine the shape and scale of the selenodetie network, 
the stellar photography and the state vectors have to be freed from their depen- 
dence on their given system of reference. This is possible with the use of minimal 
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constraints particular to each type of data. For example, instead of absolute 
camera orientation, one can use relative camera orientation, i.e., the orien- 
tation of the camera axes at each exposure with respect to that at the time of 
the first exposure in each strip. These relative orientations are coordinate 
system independent. For the use of orbital support, the approach is to first' 
establish an arbitrary system for the network through minimal constraints. 

Then positional and orientation parameters with respect to this system are intro- 
duced for each pass as additional parameters. From these solved parameters 
one may detect either the random or systematic pattern in the orientation and 
position of the orbit systems . If the differences between these position and orien- 
tation parameters are statistically significant, then the approach of separating the 
coordinate systems for each arc is justified. 

6. 2 Method for the Statistical Comparison of Data 

6.2. 1 Statistical Test for the Recovery of Inconsistencies 
Between Two Data Groups 

It has already been shown that when combining lunar terrain photography 
with any of the other observational groups (stellar photography, altimetry and 
orbital support), there exists a set of parameters that can be determined (esti- 
mated). At least the shape can always be determined, in which case the estimable 
parameters are the coordinates with respect to a system established by minimal 
constraints in which the scale is also constrained. 

To generalize, consider two groups of observations G* and G 2 and a set of 
estimable parameters X, common to both groups. If X x and X 2 are unbiased 
linear estimates of X, using groups G x and G 2 , respectively, and if the groups 
are divided into 

Gi= {l 0 ,Li} , 
g 2 = {l 0 ,l 2 } , 

where L 0 are observations common to both groups, then 
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Xi - Q 01 L 0 + Q 1 L 1 » 
X 2 = Q 02 L 0 + Q 2 L 2 , 


( 6 . 1 ) 


where Q 01 > Q 02 > Qi > Q 2 are known matrices (available from the estimation) 
of appropriate dimensions . 


Next consider the matrix 

X 0 * Xj - x 2 


[I 


-I] 


Xi 

Xa 


The observations L 0 , L x , L 2 are taken to be independent with their respective 
variance-covariance matrices So t Si > S 2 • 

Then 

IQ01 Qi 


x D — [i | -q 


Qo 2 


Qi 


Lo 
L 1 
L 2 


~ [ (Q01 + Q02) 


Q 1 1 


- Qs] 


L 0 

Li 

L s 


The variance- covariance matrix of X D is 
S D = QS U Q t 


QL . (6.2) 



M 

0 

O 

O 

1 


(Q 01 + Qos ) T 

[ ( Q01 + Qos ) Qi ~ Qs 1 

0 Si 0 


Qi 


1 

0 

0 

M 

1“ 


- Q2 


= (Q01 - Q 02 )So(Qoi - Qo 2) T + Q1S1Q1 + QsSsJ Qa ( 6 . 3 ) 
Further assume that L 0 , h 1 , L 2 are normally distributed. This is a reasonable 
assumption since L 0 » , L 2 are either observations corrupted by Gaussian 

noise or their linear transformations (synthetic observations). To indicate 
this, one can write 

L 1 ~-N(fi 1 , Si) , Lj ~ N (fXs > Ss ) » L 0 ~ N (jJ . 0 > So ) 
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where 


E (L 0 ) , Mi 


E (Lx ) » Ms - E (L s ) 


Mo = 

So = E-^(L 0 - Mo) (Eo “ Mo) T j" 

S'! = e{(L x - Ml) (Li - Mi) t } 

S 3 = E (L s - (is) (L 3 - Ms ) T j" 

E 0 is the expectation operator and £ 0 , £ x , £2 are assumed to he the true 
variance-covariance matrices of L 0 , L x , L 2 . 

Since X D is a linear transformation of L 0 , L x , L 2 , it is also normally 
distributed, thus 

Xp ~ N (Mo * Sp) 

Since Xi , X 2 are unbiased estimates of the same parameter vector, 

e{x x } = e{x 2 } , 

and therefore, 

Md = e{x 0 T = e{x x - X 2 } = EjXi} - e{x 2 } = "0 

and 

X 0 ~ N(*0 , E d ) . 

At this point a statistical test is performed on the hypothesis: 

H:X D is a sample drawn from a population with the distribution N(0 , £„). 

The statistical test is performed at a certain significance level, i.e. , the hypothesis, 
that X D belongs to the population with distribution N{0 , £ D ) with a probability 
a, where O^a^l, is tested. The hypothesis is certainly rejected for a - 1 
and certainly accepted for a. = 0. The choice of ol, from the interval 0 < a < 1, 
is subjective, but it has to be close to 1 to be meaningful. Standard choices 
are a = 95/100 and a = 99/100. For the performance of the test see [1], 
especially chapter 5 . 

Since X D is intended to be an estimate of a null vector of quantites, it is 
known a priori that the above hypothesis H is ture, and its rejection in the perfor- 
mance of the statistical test would correspond to a rejection of the estimate X 0 
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In the latter case one would have to 


or its variance-covariance matrix L 0 . 
accept the fact that the observations L 0 , Li , L 2 and their statistics £ 0 » 

Si , S 2 are statistically inconsistent at the significance level of the test. The 
reason for this inconsistency can be either a systematic error in at least one 
of the observation sets L 0 , , L 2 or one or more erroneous variance-co- 

variance matrices from £ 0 > Si » £2 • Systematic error means that an actual 
observation does not correspond to its expected value Included in the model. 

6.2.2 Testing for Inconsistencies due to Systematic 
Errors in Altimetry and Stellar Photography 


For die particular problem in this paper, consider die following set of 
observations: 

G t = jlunar photographic observations j- 
[ stellar camera observations! 


Ga 

G 0 


= { s 

= jaltimeter observations j- 

= | a priori estimates of camera exposure station 

coordinates from orbit analysis! 


Assume that group G L is free of systematic errors and that the corresponding 
variance-covariance matrix does not deviate significantly from its true value. 
This assumption is reasonable in view of previous photogramme trie experience. 

Next consider the possibility of systematic errors in groups G s and G A . 
Testing the groups, Gi = G u and G 2 = G L U G s , one can detect systematic 
errors in G s . The common parameters used in the test will be the network 
coordinates with respect to an arbitrary system established through minimal 
constraints specifying position and scale. 

In a similar way one can test the groups Gi = G u and G 2 = G L U G A . 
However, from these tests one cannot detect significant systematic errors, 
such as a scale factor in altimeter observations or a constant error in orien- 
tation of the stellar camera. To detect such errors one needs comparison with 
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orbital support, and for this purpose one should test groups like G x = G L U G s , 

G 2 = G l Uq 0 for orientation; or G a = G U U G A , G 2 = G U U g 0 for scale. 

Although such tests can detect orientation or scale errors, they cannot distin- 
guish between stellar camera errors and altimeter observations , or errors in 
the orbital support. This is a definite limitation of the system. 

6.2.3 Testing for Inconsistencies in System Orientation 
and Positioning among Different Orbit Passes 

As mentioned earlier in view of the geometry of earth-based tracking 
used in orbital analysis, it is expected that shape and scale for each arc is 
well determined, but there is a possibility of weak determination in positioning 
and especially in orientation. To investigate this matter one can test the groups 
Gi = G l U Gq UGq and G 2 = G L U GqU Gq , where and Go are given orbit 
coordinates for the i and j th arcs; while Gq and Gq are coordinates of the 
arcs where -the possibility of a change of- coordinate system has been introduced 
in the model. If no inconsistencies are justified (i.e., our hypothesis is not 
rejected), one can proceed in testing the possibility of different scales in altimetry 
and orbital support. The appropriate test for this is between the groups 

G i = g l U ( U Gq) and G 2 = G t U ( U Gq)Ug a . 

1-1 i = i 

In case the arcs don't fit together one would rather test 

G i = G l U ( u Gq) and G 2 = G L U ( u Gq)Ug a . 

i= i i^i 

Finally, if there is no orientation inconsistency between different orbit arcs, 
it is possible ;to test for stellar camera orientation only, in which case one tests 
the groups 

n n ‘ 

G 1 = G l U( U Go) and G 2 = G L U ( U gJ)Ug. . 

1 ~ 1 1=1 

At this point it should be obvious that in the search for statistical inconsistencies 
between the various types of data, one is faced with an enormous, although not im- 
possible, computational task. The problem obviously calls for a simplified, not 
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completely rigorous approach. The main problem is to determine whether 
or not the various orbit arcs are consistent in position and orientation. After 
this question has been answered, one has to find a way to take care of the incon- 
sistencies between the two groups, namely the terrain photography with its 
supporting observations <G L U G S UG A ) and the estimates of orbit coordinates 
(G 0 ) on the various passes which may be constrained to the same reference 
system or not. 

Before continuing with this problem one should note that if inconsistencies 
in system definition between the different passes are found, one must abandon 
the hope of positi oning or orienting (except via stellar camera observations) 
the network in a meaningful way. The simplest way to check reference system 
consistency between orbit arcs is given below. 

First of all, a solution must be made using terrain photography and all 
passes to determine coordinates with respect to the system specified by one of 
the arcs . All other arcs are considered to refer to a different system Which 
differs from the first one by a vector 6 1 , where 6 \ - [ 6 j T s | 6 /p] corresponds 
to the i th arc (i = 2, 3, • • • n) and 6 ls , 6 j P are vectors of shifts and orientation 
angles included as parameters and solved for. 

Next, a series of tests are performed on the hypotheses as follows: 

H ls : 6is belongs to the population N ( 0 , L ) 

0 is 

H 1R : 6 1r belongs to the population N (0 , £ ) 

0 1R 

for i = 2, 3, • * -n 

where Ss , Sc are the variance-covariance matrices of 6 ls and 6 lfl , also 
°is ir , 

computed during the .solution. 

If all hypotheses are accepted then one can obviously constrain all the arcs 
to the same system. If the majority of the tests is accepted than the coordinate 
system of these arcs is used as the reference system and all the remaining arcs 
are free of the coordinate system. 
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6 . 2.4 T esting for Ineons istenc ies betwe en Photogrammetric 
Control and Orbital Support Data 

After the problem of the different passes has been settled, one can return 
to the problem of taking care of inconsistencies between photogrammetric control 
and orbital support. Returning to the test of the hypothesis X D ~ N ( 0 , S D ), 
note that for every significant level a < 1 , one can find a scalar a such that the 
hypothesis X D ~.N ( 0 , Sp), where Sp = aSp, is accepted. 

It is obvious that an increase in S£ makes the vector X D a more possible 
candidate as a member of the population N ( 0 , S £, ) . • 

If the test refers to two groups of observations of. the form G x = L x and 
G 2 = -j^Li , L2 j- , then 

^ 1 _ Q 11 L 1 

X2 = Q21 1*1 + Q2 L2 

and Xp = X x - X2 = (Qn - Q2 x)I j x + Q2 S2 ” [ (Qu - Q2i)|Q2l 

Setting (Q n - Q 2 i) = Qi , one gets 
Xp “ QiLx + Q2L2 
and So ~ Q1S1Q1 Q2 ^2 Q 

One can always find two constants cri , a 2 such that the hypothesis ' X D ~ N ( 0 , S £ ) 
is accepted where 

Ed ” Q x (o’! Si) Qi + Q2 <cr 2 S ) Q 2 

This indicates that the effect, of the inconsistency in the data can be reduced by 
updating the variance-covariance matrices . Obviously this cannot be done in an 
arbitrary way and a proper way has to be found, conforming with some criterion 
of optimality. 

In simple terms, one needs to find a proper way of weighting orbital support 
versus lunar photography for a statistically consistent solution. This problem 
is investigated in detail in the next chapter. 
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6.2.5 The Use of Minimal and Inner Constraints 


The previous discussion repeatedly mentioned the use of a set of minimal 
constraints as a means to find station coordinates even though the observations 
may be position, orientation or even' scale invariant. 

Minimal constraints are a very useful tool, but there are some computational 
problems associated with their use. To understand this, consider an analytic 
phototriangulation for a strip of photographs where the camera coordinates and 
orientation angles of the first camera exposure point in the strip are constrained 
to a set of constants, and scale is provided by some other means, e.g., altimeter 
observations. It is obvious that in the solution with coordinates as unknowns, 
the variances of the first camera point will be zero. Furthermore, the variances 
of point coordinates will grow larger for points farther and further away from the 
first point. Since the variance-covariance matrix of the unknowns is the inverse 
of the coefficient matrix in the normal equation's, any difficulties in inverting 
such a matrix will be present in the inversion of its inverse, i.e. , of the normal 
equations' coefficient matrix. A matrix that has diagonal elements which grow 
systematically along the. diagonal, is bound to be difficult to invert, compared to 
a matrix -that has a smaller possible variation along the diagonal. 

Since a different set of minimal constraints will result in a different coefficient 
matrix in the normal equations, it is obvious that one should use a set of minimal 
constraints which is optimal from a computational. point of view. When the optimality 
criterion is the minimum trace of the variance-covariance matrix (inverse of the 
normal equations' coefficient matrix), the corresponding set of optimum minimal 
constraints is the so-called "inner constraints . " This concept of inner constraints 
has been investigated in [2], with emphasis on range observations. 

Appendix A is an extension of .this concept for photogrammetrie application. 

The use of an inner set of constraints instead of a more easily identified set of 
minimal constraints causes no additional difficulties . Formulas for the transformation 
of solution vectors and variance -covariance matrices when different sets of minimal 
constraints are used, are found in [14] . 
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7 . OPTIMAL WEIGHTING OF ORBITAL SUPPORT DATA 
VERSUS LUNAR TERRAIN PHOTOGRAPHY 


In this chapter a way is sought to use two sets of observations together with 
their variance-covariance matrices, when the two sets are statistically inconsis- 
tent (at a certain significance level). The only rigorous answer, of course, is to 
reject either or both sets of observations, the former if there is reason to believe 
that the inconsistency comes from that set only. Such a negative attitude does not 
solve the problem, and instead one should look for a way of using the observations 
despite their inconsistency. 

If there is no way of improving one or both sets of observations, then the 
effect of the inconsistency can be reduced by updating their variance-covariance 
matrices. The simplest way of updating thses matrices would be to assume that 
they differ from their consistent counterparts by scalar multiplications only. This 
way the relative accuracies between the observations within each set would be pre- 
served, while the two sets are properly weighted against each other by giving the 
less consistent set a smaller weight, thus reducing its effect in the solution'. 

More explicitly, given two sets of observations Li and L 2 » with the corres- 
ponding variance -covariance matrices K x and K 2 , it is assumed that their true 
variance -covariance matrices are 'L 1 = ajKj and Z 2 = a 2 K 2 , where the scalars 
CTi and CT 2 are to be determined. Since there are an infinite number of such 
scalars which would reduce the effect of inconsistency in the two sets of observations, 
some criterion is needed to determine the optimum pair of values for- 0 ^ and cr 2 . 

In the following, two such criteria are discussed, namely the well known Maximum 
Likelihood and the less familiar but in this case computationally more practical 
Minimum Norm (M3NQUE). 
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7. 1 Maximum Likelihood Criterion 


Lemma 1 [1] 


If x = (X! , x 2 » ' ’ 'x n ) 7 is a random vector variable with expectation E(x ) = fx , 
and its variance-covariance matrix is 

e|(x - fx) (x - u) T | = V , 

then assuming that x is multivariate normally distributed, the j oint probability 
density function is given by 


f(x 1} x 2 , 


) = 


(2ir?Vjvl 


exp 


j- J 


( 7 . 1 ) 


7.1.1 The Regress ion Model 


Suppose that a set of parameters X a and two sets of observations and L a2 
are related by the vector model f (X 8 , L aI , L a2 ) = 0. After a Taylor series 
expansion and neglection of second or higher order, the following set of linear 
equations is obtained: 


0 - f (Xq> L^i, Lbs)"*" 


5 f 
BX a 


(X a -X 0 ) + -rr 

5 L a l 

0 

b 


(L al ” L bl ) + 

0 

b 


+ 


df 

9 L a 2 


(L a 2 “ L b2 ) > 

o 

b 


where X 0 are approximate values of X a ; L w , L b2 are observed values of 
L^ , L^ and the partial s are evaluated for approximate and observed values. 
The above equation maybe written as 


0 = -W+ AX + Bx V x + B 2 V 2 

W = AX + B a Vi + B 2 V 2 . (7-2) 
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The following information is also known: 


E(Vx) = 0 


E(V 2 ) = 0 

Var (Vx ) « 

E -f V x V x T 1 

- O-fKx 

Var (V 2 ) = 

e{v 2 v 2 t } 

= <4k 2 


where Var (y) stands for the variance -covariance matrix of vector y; E-j j- is the 
expectation operator; K x , K 2 are known positive-definite matrices and a x , a 2 
are the unknown scalars to be determined. 

Some of the obvious results are 

E (W) = E(AX) + E (BVi) + E(BV 2 ) = AX . 

Also: 

Var (W) = E { (W - A X) (W - AX) T \ = 

= E ' (BxVi + B 2 V 2 ) (B 1 V 1 + B 2 V 2 ) t '1 - 

J 

= B x e{ViVI)bI .+ BiEj V!V 2 1 bI + 

+ B 2 E<V2Vx1b 1 + B 2 E v 2 v 2 t b 2 


If V x and V 2 are uncorrelated, we have 

Var (W) = BxVar(Vi)B][ + B 2 Var(V 2 )B 2 = 

= O’ 1 B iKx bJ + O 2 B 2 K 2 B T 2 = 

2 2 

- CTi Sx + CT 2 S 2 = 

' 2 2 

= ax (Sx + y 2 S 2 ) = aiH 

/ ff 2 \ 2 

where y 2 = • 

2 2 A 
If ax and cr 2 were known, then the least squares estimator X , of 

X, is given from the solution of the normal equations 
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[A T (a?S 1 +.CT|S 2 )‘ 1 A] X = • A T (criSi + ai Ss )* W 


A 

• and X minimizes the quadratic form 

(D - ~~t- VjKi V x + Vz K2 1 V 2 


7.1.2 The Maximum Likelihood Solution 


The problem is now to find from all possible values of o? and or| , the optimum 
pair, where the optimization criterion remains to be established. Recalling the 
definition of the vector W, 

W = -f (Xo»- Lb) = [ Lbi 1 ] 1 


and noticing that W is a deterministic function of the random observations L b , 
the established criterion of optimization is the "maximum likelihood" criterion, 
where the likelihood is referred to the vector W. 

The likelihood of W is given by L = f(w x , w 2 , * * *w „ ), where f is the joint 
probability density function of the components w 1 , w 2 , • • -w n of the vector W. 
Equation (7.1), after replacing x with W, fx with AX and - V with cr? H, yields 

L = exp [~ (W - AX) ’ H £ (W - AX) -| 

(2ir jVl.trf HI L 2 j 

To get a maximum of L we can alternatively get a maximum for 

X = log L , 


md X = -log<2^ - logaj - logyfif - 

4 0* x 

f 

where V| oiH | was replaced by cr£ V |H j. 

The function X obtains a maximum for these values of X, Oi , y 2 which 
satisfy the following equations 


3X_ = n 3 X 


0 


ax 


0 . 
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The developed form of these partial derivatives are as follows: 


5 A 

ax 


^[(W-AX) t H' 1 (W- AX)] = 


“ T “T [W t H _1 W- 2W t H' 1 AX+ X t A t H' 1 AX] 

ax 

1 


F [ - 2 W T H' 1 A + 2X t A t H’ 1 A] = 0, 


2oi 


T A T Tj -1 


where the following identities have been used: 

S (Ry) = R and (y T Ry) = 2y T R 


ay 


ay 


ax _ 


For cr i 0 and Oj < 00 , — = 0 gives 

o X 


X = (A T H X A)’ 1 A T H' 1 W 
The derivative with respect to cri is 

IF = T \ -^ !j -[(W-AX)'H- 1 (W-AH)] = 


n 


cri CTi 
From where for or* £ 0, one gets 


+ — §— (W - AX) 7 "H (W - AX) = 0 


2 (W - AX) T H (W - AX) 
cr 1 = 

n 


At this point, before proceeding further, two algebraic lemmata are introduced. 


Lemma 2 


For the square matrix H = H(y), the following relationship holds: 


B [ log j H | ] = tr ( H 1 ~ ) 


a y 


ay 


(7.3) 


(7.4) 


(7.5) 
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Proof: 

Based on [5], page 266, eq. (10. 8.22), for a square matrix H 
a |H ! . 3H T 


tr [ H* ■ 

3 y 1 by 

* 


] 


where H * is the adjugate matrix of H (i , e. , H jj = cofactor of H ls ) . From [19], 
page 39 

H H* = | H 1 1 


where I is the identity matrix, or H* = |h|H'\ Therefore, 

a ajHj _ i . — an T 


ir [l 0 K i “ 1 1 = wir = i 


HI trlH * ay -1 " 


1 tr[|H|H 1 4rH = tr[H 4 ] . 


Th 


ay 


3y 


Lemma 3 


For the square non-singular matrix H = H{ y ), the following relationship 
holds: 



3 y 

- -H' 1 -—- H' 1 . 

3 y 

(7.6) 

Proof: 

K- Hl = I . ” 

= o^ H H , + H aiHli 

3y 3y 



H mix 

by 

aH 4 , ^(H- 1 ) 

b>y H d 3y 

— II *1 ^ ^ JJ ”1 

ay 


□ 

The partial derivative of X with respect to y 2 is 

*ayT = ~ T ~by~ 2 [1 ° g I H 1 1 " Icrf (W ' AX > T * (W - AX) = 


= - — - trace (H 


3H 

ay 2 


) + 


3 H 


(W- AX) 1 H* (W - AX) 
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Since H ‘ = S x + y 2 S 2 » 

_^L = o 

-N fa 2 

By 2 and 

^A_ = _ -|~trace(H- 1 S 2 )+-^ 2 (W - AX) T - AX) = 0 

or, (W - AX) t H 1 S 2 H " x (W - AX) = cr 2 trace (H ' X S 2 ) . (7.7) 

If in the above equation X(y 2 ) and cr x (y 2 ) are replaced by their respective 
values from equations (7.3) and (7.4), an equation is obtained withy 2 as 
the only unknown. 

If it is known a prior ily which observations set is less reliable, one can 
always set 0 < y 2 < 1, and an approximate solution can be obtained by iteration 
on equation (7.7). 


7.1.3 Solution with Only one Updated Variance-covariance Matrix 

Consider a least squares adjustment of the linear model W = AX + BV 
where E(V) = 0, Var(V) = £, the well known solution of which is: 

X = (A t M^A)’ 1 A t M’ 1 W 

VarX = ~~ (A t M~ 1 A)' 1 

where f is ihe degrees of freedom and M = B S' 1 B T . 

^ „ - 

Consider the same adjustment, but with Var (V) = E - k S , where k is a 
scalar, then 

vs 

X = (A T M ‘ 1 A ) ’ x A t M 1 W , 

_ 1 

but M = BS B t = k 2 BSB 1 = k 2 M, M' 1 = M 1 . 
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Consequently 


X = k 2 (A t M ' 1 A ) * ( ~pr ) A T M W - X 


and also since V is a linear transformation of X, V= V. Also 

r T v " l " tr - - — tr T -3- 


=. Vfl it-rA!. , ? T k 2 fA T M* A) 

f 


Var(X) - f 

= Var(X) 


i.e., a scalar multiplication of the variance-covariance matrix of the 
observations has no effect on the results. 

Returning to the model with two sets of observations 


W = AX+ BiVi + B 2 V 2 = AX + |;b 1 B 2 ] 


Vi 

v 2 


= AX + BY 


and 


Var(V) = 


Var (V x ) 

0 


CTiKi 0 

0 

Var (V 2 ) 


0 ct 2 K 2 


Instead of Var(V), the following variance -covariance matrix can be used without 
effecting the results : 


k Var (V) = 


k 2 a?K! 


0 


k 2 a|K 2 


Selecting k = — s - * one S ets 
0*1 



Ki 

0 


Ki 

0 

Var (V) = 

0 

c 4 K 

a? *• 


0 

yK z 
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This is exactly the same as using 


Var (Vi) = K x 
Var (V 2 ) = yK 2 


where only the constant y remains to be determined . In this case 
Var (W) = BxKiBj + y Ba K 2 B<!> = Si+yS 2 = H . 


The likelihood of W is given by the same expression as before . 


X = logL = - log (2 it) - log J |H 


(W- AX) T H ' 1 (W- AX) 
2 


Setting ~~r = 0 and 

o X, 


= 0 , one can derive that 


(A T H' l A)X = A T H*W 


(7.9) 


(W - AX) T H ' 1 H ' L (W - AX) = trace (H’ 1 -™-) 

ay dy 


Since H = S x + y S 2 > 


— — = S 2 , so that 
3y 


(W- AX) t H' 1 S 2 H'^W-AX) = trace (H 1 S 2 ). 


(7.10) 


H ' X S 2 H * = (HS 2 H)' 1 = a* 


O = HS|H = (S x + yS 2 )S'i (Sj + ys 2 ) = 


S 3 .S 2 Si + 2 y Si + y S 2 


Also 


H’ 1 S 2 ~ [S 2 1 (Si + y S 2)] 1 = (S 2 1 Si + yl ) . Substituting the above 


equations into (7.10), one arrives at the equation which is to be solved for y 


(W - AX) T O ' 1 (W - AX) = trace [(S 2 1 Si + yl ) * ] 


(7.11) 


where 


X = (A T H A ) A T H W 
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and 

a = SiS^Sa + 2y S i + y S 2 . 

If y is needed with alow accuracy only, and the algorithm is arranged so that 
0<y< 1, equation (7.8) will provide X for a given y . The value of y at which 
X is maximum will also make the following X' function a maximum: 

A ' - -log/liT- tw-AX)^ <w . -.«L 

Thus computing X' for various values of y, at some equal increments (say, 
for example, y= 0. Ik, with k - 1, 2, • • -10), a graph of X' = X'(y ) can be 
constructed, which at its maximum will provide the sought y . 

7.2 The MINQUE Criterion 


The criterion for optimal weighting in this case is the minimization of an 
appropriate norm, established by C. R. Rao in his so called MINQUE theory 
(Minimum Norm Quadratic Unbiased Estimator). The theory is developed in 
[15], chapter 4j, and in a series of papers ([16], [17] and [18]). Only a very 
short outline of the theory is given here. Rao considers the linear model 

Y = X/3 + c 

where, e = Ux^! + U2I2 + U k 4 k = U£ in which 

E {^ 1 } ~ ° ^d = o - i6 1J I (Cj are unknown) 

E(e e T ) = (TiUjUl + •• • CTkU k Uk = o\ V 1 + • • • ex* V k , 

where Y are the known observations, X is a known matrix and /3 is a vector of 
parameters to be estimated. The quantity to be estimated is a linear combination 
of the CTi 

SPiffi = P T cr 

i 

where p T = [p 1} P 2 > * ’ ■ Pk 3 and a T = [erf , ai » • * * o* ] . The proposed estimator 
is a quadratic form Y T AY, where the matrix A remains to be determined. The 
estimator is also required to have some additional properties which impose 
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additional constraints on the matrix A. These properties are the following. 

a. ) Invariance . If j3 ' - j3 - j8 0 , where /3 0 is a constant vector, the following 
relationship needs to be fulfilled: 

Y /t AY' = Y t AY 

where Y ' - XjS ' + e . This property imposes the condition that AX = 0 . 

b. ) Unbiased . The following relationship needs to be fulfilled: 

e{y t Ay}- = SpjCTi = p 7 a . 

The corresponding condition is that Trace(AV t ) = p x . 

c. ) Minimum Norm . If the £ t 's were known, a natural estimator would be 

a i 

Jli 

where n t = number of elements in £ t . Then £ p j o \ = E £j£i = £ T A£, 

A i i n i ' 

where £ T = [£ { ,. £2 . * * * £k] and =6 ^ . 

n j 

The difference of the estimators is 

Y t AY - £ t A£ 

from where, using the condition AX = 0, one gets 

Y t AY- £ t A£ = £ t U t AU£ - £ t A£ = £ T (U T AU- A)£ s 
It is desired to have the norm ||u T AU - A || = min. for some properly defined norm. 

7.2.1 Solution Under Euclidean Norm 

Consider the Euclidean norm defined as || B|| 2 = SS(B tJ ) 2 , then the norm 

i i 

to be minimized will become 

||u T AU - A|1 = Tr(AVAV) - Tr(AA) 

where k 

V = S Vi = S UiUi T . 

1=1 i 

Since A is fixed, the problem reduces to minimizing Tr(AVAV), subject to the 
conditions AX = 0 and Tr(AVi) = p t . The solution given by Rao is as follows: 
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A — S X[ RV i R 
1 

where R = V 1 - V'^X 1 V^X^X 7 V ' 1 and X t is the solution of 

SXiT^RViRVj) = Pj 

i 

or S X - p 5 

where X T - [Xi , X 2 » * " X k ] and S 13 = Tr (RV j RV 3 ) . 

For the estimated quantity 

p T q = SPia? = Y t AY = X t Q, 
i 

-where Q = [Q x | Q 2 | ••• Q k ] and Q t = Y T RViRY. Combining SX = p and 
p T a = X t Q, one gets 

& = S -1 Q 

where again S u = Tr(RVj RVj ) and Q 4 = Y T RV t RY. 

7.2.2 Application of the MINQUE Theory to 
the Optimal Weighting Problem 

The model is (see equations following (7.2) ): 

W = AX + BiVj + B 2 V 2 = AX + e 

E(V 1 ) = E(V 2 ) = 0 E ( V x V 2 ) = E(V 2 Vl) = 0 

E (V x V 1 ) = a?K x E(V 2 VJ) = af K 2 

Var(W) = E (ee T ) = B a E(ViVl)Bl + B 2 E(V 2 V 2 )B 2 = 

= erf BxKxBf + erf B 2 K 2 b£ = 

~ ax Ex + C 2 E 2 . 

The MINQUE solution for af and cr 2 is 
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where 


S« = Tr(RS, RSj) , q 


with 


Qi 

q 2 


q i = W T RS t RW , i = 1 , 2 
R = H' 1 - H' 1 A(A t H 1 A)' 1 A t H ' 1 
H = Si + S 2 . 

The computational algorithm is summarized below: 

Sj = B3.K1B1 

S2 ~ b 2 k 2 b| 


H = S x + S 2 

R = H ' 1 [I - A(A t H' 1 A)' 1 A t H' 1 ] 
q 1 = W t RSiRW 
q 2 = W T R£ 2 RW 

511 = Trace (RSxRSi) 

512 = Trace ( RSi RS 2 ) 

Sai = Trace (RS 2 RS x ) — 

S22 = Trace (RS 2 RE 2 ) 


and 


af 


Su 

Sj.2 

-1 

qi 



S 21 

s 22 


q 2 

_ _ 



- 


^ _ 


D = CT1E1 + cr| S 2 


and 


X = (A t D' 1 A)‘ 1 A t D‘ 1 W. 


The main computational task is in the computation of R, but the computations 
are of the same order as those involved in a normal least squares solution 
when oi and a I are known. 
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8. SUMMARY AND CONCLUSIONS 


This work consists of two problems in the establishment of selenodetic controls 
the determination of the shape and scale of a network of points on the lunar surface 
and the determination of the relationship between this network and a physically 
meaningful selenodetic system, namely that defined by the selenocentric principal 
axes of inertia. The main source of shape information is lunar terrain photography. 
Supporting altimetry and stellar photography, beyond the obvious scaling and orien- 
tation, also contribute significantly to the strength of the network. The role of 
coordinate systems, especially those inherent in the external information (lunar 
theory, physical libration theory), has been analyzed and it has been shown that 
positioning and orientation of the network is possible only with respect to a.) the 
estimated "selenocenter" provided by the lunar theory used in the orbital analysis 
and b.) the estimated "principal axes of inertia, " provided by the libration theory 
used in both orbital analysis and the reduction of stellar photography. The possi- 
bility of orientation of the network with respect to the principal axes of inertia 
system by means of gravity information has been shown to be impossible in view 
of the present limited accuracy of second degree harmonics of the gravity field of 
the moon. In view of indications in previous work of inconsistencies in the esti- 
mates and statistics of altimetry, stellar photography and orbital support data, 
statistical tests have been described for the recovery of such inconsistencies. 

The role of minimal constrained soltuions has been emphasized in retaining 
coordinates as parameters, although the observations involved in solutions and 
statistical tests may not provide information on positioning or orientation or 
scale. 

An algorithm has been developed for the use of inner constraints, as a compu- 
tationally optimum set of minimal constraints . With special reference to the orbital 
support, the possibility of inconsistencies in system definition between different 
passes, has been elaborated. Procedures for statistically establishing such 
inconsistencies have been outlined and the idea of the utilization of such passes 
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as free of positioning and orientation has been introduced as a means of reducing 
the effects of their inconsistencies. 

Finally, in the belief that inconsistencies in orbital support data are partly 
due to non-realistic statistics, arising from the uncertainties in lunar gravity 
(on which orbital analysis heavily depends), methods have been introduced for the 
optimal weighting of orbital support versus terrain photography. No computational 
effort using real data has been made, however, towards the use of the outlined 
techniques. It is therefore proposed that the following numerical investigations 
be undertaken: 

1. ) Lunar photography with orbital support data should be used in a 

solution with respect to an arbitrary system, established by minimal 
constraints. In this solution each orbit pass should be treated free of 
position and orientation, and the estimates of shifts and orientation angles 
of each orbit pass and their variance-covariance matrices should be tested 
for statistical significance. 

2. ) A second solution should involve lunar photography and orbital 

support data, where only those passes with significant inconsistencies 
are treated as coordinate system free. This solution should be com- 
pared with solutions in which altimetry and stellar photography data 
are included, and the results should be tested for inconsistencies in 
those latest data sets. 

3. ) If no inconsistencies have been found in altimetry and stellar 

photography, or if they have been removed, a final solution would 
involve all data types, where orbital support is optimally weighted 
versus the lunar photography data, for the removal of any remaining 
inconsistencies . 
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Consider the photogramme trie problem: L = f(X x , X 2 )> where L are 
photo coordinates of a cluster of points; X x are the unknown coordinates of 
those points and X 2 are the unknown elements of the exterior (and possibly 
interior) orientation of the camera. 

The application of a least squares adjustment leads to a set of normal 
equations 

NX + U = O 

where the matrix N is singular with rank deficiency seven, due to the lack 
of orientation positioning and scale in our model. If N is a u xu matrix, then 
.one way of obtaining a solution for X is to introduce a set of minimal constraints 
70^!= 7O ! , such that 

N E t 

det [c D 

Among all possible C matrices there exists some matrix E such that if 



then Q is the pseudo-inverse of N, i. e. , a matrix fulfilling the four relations: 

NQN = N 
QNQ — Q 
NQ = (NQ) T 
QN = (QN) T . 

The set of constraints EX = 0 is then said to be a set of inner constraints. 

It can be shown that the matrix E also fulfills the relations 

det (EE T ) ^ O and AE T = 0 
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or, since N = A T PA, the equivalent NE T = O. 

The first step in the investigation to determine the matrix E was to follow the 
procedure in [2]. According to this procedure, the matrix E is formed as the matrix 
of the partials of the parameters X, with respect to differential translations da x » 
da Y > da z , differential rotations d<p x , d(p 2 , d<p 3 and differential scaling d e 
of the coordinate system. The results showed that the changes in the angles of 
camera orientation were non-linear functions of the differential rotation components . 
Following the approximations in [2], page 18, 


cosd^jRi 1, sin dcp ! « d(p t , d<p x - d<pj = 0 


a linear relation is obtained but the rsulting E matrix does not iuiiui tne relation 

7 

AE = 0 . A more careful investigation showed that the above approximations are 
not valid and that Blaha’s results, which are correct despite the approximations, 
can be derived rigorously. 

The question that posed itself next, was whether or not a set of constraints 
that does not involve the camera orientation parameters can be an inner set of 
constraints . Augmenting the normal equations 


p 

- 


— ■“* 


r “ 

Nn 

Nia 


X x 

+ 

u x 

N i2 T 

N 33 


X 2 


u 3 


- 


- 


- - 


in such a way that X x are the coordinates of the control network and X 2 the 
parameters of camera orientation, a set of constraints involving only X x can 
be written as follows: _ 


X x 


Ei Xj = 0 or EX = [ Ei 0] 


X ? 


Recognizing that N - A PA or for simplicity and without loss of generality 
N = A T A > one gets 


N - A'A = 


Ax 


[Aa A 8 ] = 


AiAj. A^Ajg 

aJa s aJ^-2 


Nix N ls 
Nl2 N 22 


Now select a matrix E x such that 
det (E X E^) and 


T T 

N XX Ex = 0 or A X E X ) = 0 
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Then the matrix E fulfills the relations 


EE = [E x 0] 


Ei T 


EjE l det (EE T ) = det (E x E\ ^ 0 


T 

NE 


Nn 

N xs 


e x t 


N n Ei T 


1 

t> 

w 


' 0 " 

nJ 

n 23 

J 


0 


fc- H 


1 

tH 

W 

<5* 

h W 

c 

1 


0 


and therefore EX=0 is an inner set of constraints. 
The constrained normal equations now become: 


- - 


“ 


- 

Nn N 1S El 


Xx 


Ui 

Ni 2 N S2 0 


X 2 

+ 

U 2 

E x 0 0 


K c 


0 

. — 


* 




where K c is the vector of the Lagrange multipliers. If the order of equations 
and unknowns is changed 


T 


* ~ 


— 

Nn Ei N ia 


X x 


Ui 

Ei 0 0 


Kc 

4- 

0 

\ 

£ 
<0 ^ 

o 

£ 

L~_ 


x 2 


U 2 


Since N is singular both N u and N 22 must be singular, because otherwise 
N could be inverted by partitioning utilizing the inverse of the non-singular 
matrix among N ai and N 22 . Now it is possible to invert by partitioning, 


Nn 

e t I 
Ei J 

0 1 

Nl 2 


r 


_E_i__ 

0 

since 

Nn El 

is not singular. 

Ni T 2 


n 22 


Ei 0 

, — 



The remaining question is how to find the matrix Ei in such a form that 
A x eJ = 0 (N n Ei = 0) and det (E x E^) ^0. 

If N X1 is an nxn matrix it can be augmented 


N X1 


Nin 

Nil 2 

n-7 x n-7 

n-7 *7 

T 


n x22 

7 x n-7 

7x7 
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♦ 

in such a way (by rearranging the order of unknowns if necessary) that N U1 
is a non-singular matrix and N 122 = N lls N U i N 112 . 

Setting also E x = [ E X1 E 12 ] one must have 

F X N F X N-F F X F 


Nu Ej - 

" 

Nm 

Nu 2 


I 

J-F 

H s 
1 


Nu 2 

Ni 22 




- 

- 


- 


Nm 

Ei 1 * 

+ Nii 2 

Ex? = 0 

„ T 

T 

+ Ni 22 

r 

Nil 2 

Eh 

Ei 2 = 0 


x 

If E X 2 = I, the first of the above equations gives 


Ei? - -Nut N ll2 

and setting this value into the second equation one obtains 


Ni gg - N 112 Nni Nn 2 

which relation is known to hold a priorily. 
Therefore , 



r 



zn 

Ei 1 

— 

-Nui N 112 


_ r 




Ei 2 


I 


A different approach will be to find the matrix Ei in analytical form 
with some sort of systematic pattern, such that AjEx = 0 and det (EjEi) ^0. 
However, such an approach will not reduce the total computational effort, 
since to obtain a solution the matrix Nm will have to be inverted ex- 
plicitly or implicitly. 


Now the solution to the system becomes 


r -1 


r 1 


- - 


- 


r- 

Xi 


Nu E I N i2 


Ui 


Ql 1 Ql 2 Q13 


Ui 

N c 

=" 

Ei 0 0 


0 

=- 

T 

Qi 2 Qgg Q23 


0 



Nil 0 N 22 


_ u 2j 


T T 

Qi 3 Q23 Q33 


_U 2 _ 


If only the parameters X x are of interest, then 


X x = -Q 11 Ui -Qi3 U 2 


Setting 


• 

-1 

- 

Nu Ei 

= S = 

$L 1 2 

3 0 J 


( 

- 

CO 

*0 
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and inverting by partitioning one obtains 




Sll St 2 


'n 12 

Q33 “ ^ 

i 

N 2g - [ nJ 0 ] 

1 

J 8 

w 


0 


Q33 ~ ( N 22 “Nj, ® Su N ls ) 1 


[~Ql3 



Si 3 


N la ] 

1 

£> 

ai 

l“ 


s kl 



, 

O 

1 


Q 3i 


Q13 Sti N is Q 33 


Q11 Qi 2 


\x Sts' 

■f 

s i 1 St 2 


* 

Q 33 [nJo] 

Sli s l2 

T 

Qi 2 Q22 


St 2 s 22 


1 

p 

P 

* 


0 


Sts Sti 


Qu “^11 + Sti N l2 Q33 N 12 Sj.! 


can be found by inverting the matrix 


Nn E : 
E, 0 


Another approach is to recognize the fact that S X1 is the pseudo inverse 

4 * 

of N lx (S 11 = N X1 ) and to use some direct algorithm to find this pseudoinverse* 
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Appendix B (Referenced in section 5. 1) 


COMPUTATION OF VARIANCE-COVARIANCE MATRICES 
OF CAMERA ORIENTATION ANGLES WITH 
RESPECT TO A MOON- FIXED SYSTEM 


If K j , , £l s are the orientation angles at j th camera exposure with 

respect to an inertial system, and , 0j , €j are the Eulerian angles re- 
lating a moon-fixed system to the inertial one, then the relationship between 
the orientation angles of the camera (x } , <p } , 40 j) with respect to the moon- 
fixed system is 

M (xj.pj.w,) - M(K J ,* J ,O j )R(0 J ,e J ,€ J ) 

where 

M = R 3 (Xj)R2(<£j) Ri(Wj) 

M = BatKjJRa^jJRitOj) 

R = R 3 (-€ J )R x (e j ) R 3 (-0 j) • 

Introducing the notation 

E = [Xj <Pi co 3 ] t 
F = [K j 4 j J2 J ] T 

G = [Cj 0j] T , 

the corresponding variance-covariance matrices are given by 

Sf = Qf Sf Qf + Qg £g Qg > 

3 E , ^ 3E 

where Q F " and Qg " ~^T • 

If m j j are the elements of the matrix M, one gets 





<Pj = tan 


l ( 


m 3 i 


J 


2 

m S 2 


2 

m-33 


tan' 1 


( \ 
V m33 J 
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and: 


0 


0 


0 


mu 

m 2 i 


3E 

3a 


m 2 i -m n 
0 0 cos<pj -taxL(p 3 -tan tp } 

0 0 - m 33 0 m 3 i 


_3_ 

3a 


m 3 i 

m 32 

m33 


There is also a need to compute 3 M/ 3a , where a is an element. of F or G. 


3 M 8M ~ 

3 a F 3 a p 

where a F is a component of F, 
Using the matrices 


and 

and a G 


3M ~ 3R 

— - M — — , 

3a G 3a G 

is a component of G. 


0 

0 

0 


0 

0 

-1 


0 

1 

0 

0 

0 

1 

to 

11 

0 

0 

0 

> P'3 “ 

-1 

0 

0 

0 

-1 

0 


1 

0 

0 


0 

0 

0 


one gets the following relations: 


M = R 3 (K)R 2 ($)Ri(0) 

R - = 

R 3 (-C)Ri(0)R 3 (~ip) 




— M = P S M 

r^-R = 
3 e 

-p 3 r 

= R 3 <K)R> R 2 (0)Ri(f2) 

30 R " 

R 3 (-e)Pi Ri(0)R 3 (~tp) 



rs/ 

mT m = MP * 

* 

"H?3 t 
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Appendix C (Referenced in section 3.2) 


ON THE DETERMINATION OF THE DIRECTION OF THE 
PRINCIPAL AXES OF INERTIA FROM 
GRAVITY FIELD INFORMATION 

1. Transformations of Moments and Products of Inertia due to Changes 
in the Coordinate System 

Assume only Cartesian systems with origin at the center of 
mass of the body in question. If 0, x, y, z and 0,x,y,z are two such systems 
with corresponding moments and products of inertia, A,B,C,D, E, F and 
A, B, C, D, E, F , then the transformation from one system to the other can 
be represented as: 


, X 


X 

y 

= R 3 (x)R 2 (<P)Ri(w) 

y 

z : 

| 

z 


The problem is to find A, B, C, D, E, F given A, B, C, D, E, F and the angles 
x , (p , w . Introduce the following auxiliary coordinate systems: 


x' 


X 


it 

X 


x y 

y' 

= Rxlco) 

y 

9 

y" 

= R 2 (cp) 

y 7 

z' 


z 


z" 


z' 

H 


„ it 



X 


X 



y 

= R 3 (x) 

y" 



z 


z” 



J 


- 




The corresponding moments and products of inertia are 
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/ (z ,a + y ,2 )dm » 

* / M 

/’[(coscoy + sincoz) 2 + (-sincoy + cos to z) 2 ] dm- - 
J H 


h 


(y 2 + z 2 ) dm = A 


/ (x' 2 +z /2 )dm - f [x 2 + (-sincoy + coscozfjdm = 
** H J M 


. 2 

sin a) 


yV + y 2 )dm + cos 2 to y> + z 2 )dm 


2coscosinco / yz dm = 


X- 


sin 2 co C + cos 2 co B - sin2coD 


/ 


<x' 2 +y /2 )dm = J [x 2 + (cos coy + sincoz) z ]dm 

H M 


cos 2 co f (x 2 + y 2 )dm +. sin 2 co A** + z 2 )dm + 
J M 

2coscosinco / yz dm 

J M 


cos 2 co C + sin 2 co B + sin2co D 


f y'z' dm = f (cos coy + sinco z) (-sin coy + cos co z)dm 

M -/ m 

sinco cos co [ I z 2 dm - / y^dm] + 


(cos 2 co - sin 2 60 ) / yz dm = 
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= sinco cosw [ 


A + B - C A - B +C 


] + cos 2 co D 


E 7 


F' 


sinco cosco (B - C) + cos 2co D 


/ x'z' dm = / x(-sincoy + cosco z)dm = 

M ' J ft 

-sinco / xy dm + cosco xz dm = -sinco F + cos co E 
J K H 

/ x 7 y 7 dm = / x(coscoy+ sincoz)dm = 

J M J M 


cosco F + sinco E 


In summary 


Similarly 


A' = A 

B 7 = cos 2 co B + sin 2 co C -sin2coD 

C 7 = sin 2 co B + cos 2 co C +sin2coD (C.la) 

/ 1 

D = “ sin 2 co (B - C) + sin 2 co D 
& 

E 7 = cos co E - sinco F 

F / = cos co F + sinco E . 

A" - • cos 2 <p A' + sin 2 <p C 7 + sin2<p E 7 
B 77 = B 7 

C " .= sin 2 <p A 7 + cos 2 <p C 7 - sin 2<p E 7 (C . lb) 

D 7/ - coscp D 7 + sincp F 7 

E /y = ^ sin2p (C 7 - A 7 ) + cos-2p E 7 

f" = -simp D 7 + cosip F 7 
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and 


A = cos 2 x A" + sin 2 xB" - sin2x F" 

B = sin 2 x A " + cos 2 xB" + sin2x F" 

C = C" (C.lc) 

D = cos x D" - sin x E " 

E - sin x D" + cos x E" 

F = ■- sin2x {A" - B") + cos2x F" . 


Combining the above three sets of equations , the following final 
transformation equations are obtained: 

A = A^ cos 2 x cos 2 <p I 

B | sin 2 x cos 2 to + cos 2 x sin 2 <p sin 2 co + ~ sin2x sinp sin2co j- 


+ 

+ 


+ 

+ 


C | sin 2 x sin 2 co + cos 2 x sin 2 <p cos 2 to sin2x sincp sin2co 
D | - sin 2 x sin2 to + cos 2 x sin 2 <n sin 2 to + sin2x simp cos 2 to j- 
E -|cos 2 x sin2cp cos to - sin2x eos<p sinco | 

F -f- cos 2 x sin2<p sinco - sin2 x cos<p cos to j (C. 2a) 

t i J 


B = 


A isin 2 x cos 2 <p j- 
+ B |cos 2 x eos 2 co + sin 2 x sin 2 <p sin 2 to - ~ sin2x simp sin2toj- 

+ C |cos 2 x sin 2 co + sin 2 x sin 2 <p cos 2 t o + “ sin2x sincp sin2to j 

+ D •£- cos 2 x sin 2 to + sin 2 x sin 2 <p sin 2 to - sin2x simp cos 2 to T 

E |sin 2 x sin 2<p cosco + sin 2x cos <p sinco 1 
F sin 2 x sin2<p sin to 4 - sin2x cos<p cos to j (C. 2a) 


= A 


| sin 2 <p j 
B | cos 2 <p sin 2 to j- 
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+ c 


|cos 2 <p cos 2 co| 
D |cos 2 <p sin2co j- 
■|-sin2<p cosco} 
+ F |sin2<p sinco j- 


+ E 


(C. 2a) 


D = 


A -j^sinx cos<p sin <p j 
+ B ^sinx sin2<p sin 2 co + ^ cosx eos<p sin2col 

+ C -f- ~ sinx sin2<p cos 2 c o - ~ cosx cos<p sin2co| 

+ D j- ~ sinx sin2<p sin 2co + cos x coscp cos 2coj 

+ E sin x cos 2<p cosco + cos x smp sinco j- 
+ F |sinx cos2<p sinco + cosx sin<p cosco 1 


(C.2b) 


E 


+ B I- 


A^-cosx cos <p sin<p]- 

fl t o .l « * 

|^— cosx sm2(p sin co+ -'smx cos <p sin2coj- 
+ C cos x sin 2^5 cos 2 co ^-sin x cos <p sin 2 coj- 
+ D -j— cosx sin2<p sin 2 to + sinx coscp cos2coj- 
+ E -|cosx cos 2<p cos co + sinx sin<p sinco} 

+ F cos x cos 2 <p sinco + sin x sin<p cos co j- 


(C. 2b) 


F = A j^cosx sinx cos z <p j- 


+ B sin2x sin 2 <p sin 2 co cos2x simp sin2co sin2x cos 2 co } 

+ C sin2x sin 2 <p cos 2 co eos2x sinip sin2co sin2x sin 2 co j" 

+ D sin2x sin 2 ?) sin2co - cos 2x sin<p cos 2co + “ sin 2 x sin 2co } 

+ E j^— sin2x sin2<p cos co + cos 2 x cos <p sinco} 

"[-~sin2x sin2<p sinco + cos2x eosip cosco} (C.2b) 


4- F 
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2. Determination of the Directions of the Principal Axes of Inertia 


Assume that the system x , y , z is the principal axes of inertia system, 
and x , y , z is an arbitrary system with the same origin for which 
the integrals A , B , C , D , E , F are known, then 

D = D (A, B, C,D,E, F, x , , co) = 0 

E = E (A, B,C,D,E,F, x , <p,CO) = 0 

F = F (A, B,C,D,E,F, x ,<p , CO) = 0 


This is a system of 3 nonlinear equations which can be solved for the 3 unknowns, 
x , <p , co . The 3 angles then, specify the orientation of the principal axes of 
inertia with respect to the arbitrary system x , y , z. However knowledge of 
all the 6 quantities A, B,C,D,E, F is not necessary. A little algebraic 
manipulation of the transformation equations for D, E and F results in 


D 


A - B 
2 

A - C 


|sinx sin2<p sin 2 co.- cos x cos<p sin 2co j 
sin x sin 2<p cos 2 co+cosx cos<p sin2coj- 
+ D ~ sinx sin2<p sin2co + cos x eos<p cos 2co j- 
+ E j^~sinx cos 2 <p cosco + cos % sin<p sin col 
+ F { sinx cos 2tp sinco + cos x sincp cos co} = 0 


(C.3) 


E - 


A - B 
2 

A - C 


|-cosx sin 2(f> sin 2 co - sinx cos<p sin2coj 
|-eosxsin2<p cos 2 co + sinx coscp sin2coj- 


+ D |—cosx sin 2<p sin2co + sinx cos <p cos 2 co? 
+ E ^cosx cos 2 <p cosco + sinx sin<p sincol 
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F -j -cos x cos 2<p sinto + sin x sin<p cos to j- = 0 (C.3) 


F = 


A “ B f . . 2 2 o ^ 

— - — -j cos 2 x sm<p sm2to - sin 2 x sin <p sin to + sin 2?tcos to j- 

A - C f ~ ... . _ ..... ^ .2 2 . - ~ . 2 1 


A - C r 2 2 2 

+ — | -cos 2x sintp sin2to - sin2x sin <p cos to + sin2x sin to j- 

+ D | ~ sin2x sin 2 ip sin2to - cos2 xsimp cos2to + sin2x sin2to| 

+ E { -^sin2x sin2'0 costo+ cos2xcos<p sintol 

+ F {-~sin2x sin2<p sinto + cos2x cos<p costo} = 0 . (C.3) 


From the above equations one can see that if the quantities A-B, A-C, 
D, E, F are known, x , <p , to can be computed. These quantities are 
related to the second degree harmonics of the gravitational field. 


From [8] (p. 160, eq. 21.043) 


A + B 


Ca-i — E 


(C.4) 


S 2i ~ D 


B - A 
4 


S 22 = “F , 


where C nt2 , S nB are coefficients in an expansion of the gravitational potential 
of the form ([8], P- 159) » 

CO n 

V 1 Y* 

v = G ) n + 1 / [C nffi cosmA + S nB sinmA] P nB (sin<p) 


r , cp , A being spherical coordinates and G the gravitational constant . 
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The inverses of the above relations are 


A 

- B 

= “ 4 C 22 

D 

H 

N 

CO 

II 

A 

- C 

_ C 20 ~ 2 C 22 

E 

" C21 



F = 2 S 32 




(C.5) 


An algorithm is still needed to solve for x»p, to in terms of C 2 o » C 2l , 

- C 22 > S 2 x , S 2? . One can always define the- x , y , z system to be close 

to the x , y , z . Then as a first approximat ion set 

cosh = cos = cos 40 = 1, 

sinx = x, sinco = <p, since = co and 

x 2 =<p 2 = co 2 = x<p = xco = <p co = ' 0 . 

Under these simplifications 

A = A + 2<pE - 2xF 

B = B-2coD+2xF 

C = C + 2coD-2<pE 

(C.6) 

D = (B - C)co + F(p -Ex + D = 0 
E = -Fco + (C - A)<p +DX+ E = 0 
F - E CO - D<p + (A - B) x + F = 0 . 

The solution of the last three equations is 


CO 


C-B - F 

E " 

'1 

D 


= 

F A-C 

-D 


E 

7t 


-E D 

B-A 


F 


And the final results are 
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U) 


DTA - a 3 1 - E F(o; + g) 
d 


(C.7) 


= E [ A + g jS - tt 2 1 + DF(2g - g) 

F[ A + a j3 - |S 2 ] + DE(2j8-a) 

H a 

where 

A = D 2 + E 2 + F 2 

a = A -B 

jS ” A - C, and 

d = (a - j3)D 2 + /SE 2 - aF 2 - o[&(oi~ /S) . 

Of course the values of x , (p , to are only first approximations but 
they can be improved by an iteration scheme based on the following 
sets of equations: 


(A-B) 1 + 1 = 


+ 


+ 


+ 

+ 


(A - B) . [cos2x - cos 2 oj t - cos 2x t sin 2 <pi sin 2 co i 
- sin 2x ; simp , sin 2co ! ] 

(A - C) j [cos 2x j sin 2 <n t - cos 2x j sin 2 <p { cos 2 to t 
+ sin 2x < sincp t sin 2co j ] 

D j [cos 2x i sin2oo j + cos 2x j sin 2 cp j sin 2co s 
+ 2 sin2x t sin<p, cos 2oo i ] 

E j [cos 2x i sin2<pj cosoii - 2sin2x! cos<Pi sincoi ] 
F t [- cos 2Xj sin 2^! sincOt - 2sin2xt coscPj coscOj] 


(A-C) i + 1 - 




(A - B)j[ - — sin2x t sincOj sin2o>i - sin 2 Xi cos 2 to i 

Li 

+ cos 2 <p t sin 2 oo £ - cos 2 x^ sin 2 <p! sin 2 to t ] 

(A - C) j [ - sin2Xi sintpi sin2toi - sin 2 Xi sin 2 toi 
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+ cos 2 <p £ cos 2 C 0 i - cos 2 Ttj sin 2 <p ! cos 2 co i ] 

+ D i [- sin 2 x { sin2(Oi + cos 2 x } sin 2 ^ sin2c0t 
+ sin 2x t simp t cos 2co j - cos 2 <p t sin2c0i ] 

+ E t [cos 2 X{ sin2<p t cos C 0 j - sin2x! cos^j sincOi 
+ sin 2<p t cosajj] 

+ F 1 [-cos 2 x 1 sin2<pi sin co t - sin2Xi cos<p! cosco- 
- sin2(D, sincoj] . 

D i + i , E 1+ i , F 1+i can be obtained from Eq. (C.3)if in the 
right hand side (A - B), (A - C), D, E, F, x, <p, co are replaced with (A-B)j , 
(A - C) t ,'Dj ,’E j, F 1} Xi , <p i» C0i , respectively. Also, X 1+ i,<p 1+ i, 
co 1+1 can be obtained from Eq. (C.7) if on the right hand sides, 

D, E, F, {A - B), (A - C) are replaced with the same quantities with a 
subscript "i". 
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3 . The Accuracy of the Directions of the Principal Axes of Inertia 

The last three of equations (C. 6 ) can he rewritten as 

D — (C20 r ^^22 ) 0 ) + 2 S 22 <p - C2 1 X + S2 1 = 0 

E - 2 S 22 ^ {2C22 “ C20 )(p S21 X 4 - C2 1 = 0 

F = C21CO - S21CD - 4C22 ^ 4 ’ 2S22 = 0 . 

Introducing the notation 

f = [D E F] t 

e = [x co to] T 

[C20 C22 C21S21 S22] 

one can write the equations determining x , <p , to as 

f = f(e, M) = 0 

and the solution will be of the form 

e = e(M). 

The respective variance-covariance matrices of M and e, E H and S e , 


s ° - [it] 

Taking the total differential of f = f(e, M) = 0 

* - [»iir] dM + [tt] de = 0 

so that 


de 

_ raf 

L Se 


r— 

Lsm 

-J dM 

f Se' 

1 = 4 

Sf 1 

1 4 r 

111 

LoM. 

J L 

SeJ 

L 

SMJ 


(C.8) 


are 


(C. 9 ) 


(C. 10 ) 
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The matrices jj^~j can be found from equations (C.3), after A-B, 

A-C, D, E and F are substituted with C 20 , C 2 2 > C 20 , S 2i , S 22 by means 
of equations (C.5). However since the angles x, <p and 00 are small, one can 
use equations (C.8) to obtain 




“C 21 


2 S 2 2 


C 20 + 2 C 2 2 

S 21 

2C 

22“ C 

20 

- 2S22 

“4 C 22 


-S21 


c 2 i 


US 

2t0 

-x 

1 

2<p 


-(P 

2(p 

1 

X 

- 2 OS 


0 

-4x 

us 

-<p 

2 



After analytical inversion 


(C.ll) 



1 

d 


2S 21 S 22 - C 2 i(2C 22 ~ C 20 ) | 2C 2 iS 22 + S 2 i(C 2 o + 2C 22 ) 

1 

C21S21" 8C22S22 I C \\ - 4 C 22 (C 20 + 2C22) 

I 

S 2 i“ 4 C 2 2 ( 2 C 2 2 _ C 20 ) J C 2 iS 2 i + 8C 22 S 22 


4SL + 4C 22 - do 
2C2lS22 ” ^ 2 1 (C 20 * 2 C 22 ) 
2 S 21 S 22 + 1 (2 C 2 2 ” C 20 ) 


(C. 12) 


where d = 2C 22 (2 c| 0 + c| x - 8 c ! 2 + S 2 21 - 8 S| 2 ) + 020 ( 8 !! - cl x ) . 

With and from (C.ll) and (C.12), one can compute 

by means of (C. 10) and finally with S M , one can find 2) e from (C. 9). 
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As an example, to get an idea of the order of the magnitude involved, 
examine the simple ease when x = <p = to = 0, C 21 = S 2 i = S 22 = 0, all 
the coefficients have the same variance s 2 without any correlation, and from 
[10], C 2 o= -204.8 X 10“ 6 , C 22 = 22. lx 10" 3 . 

For this case 

0 0 C 20 +2 C 2 2 

0 2C 22 -C 2 o 0 

- 4 C 2 2 0 0 




0 0 

0 (2C 22 -C 20 ) 

(c 20 + c 22 y 1 0 


(-4 C 22)' 1 

0 

0 


0 0 -1.13 

« 10 4 0 0.40 0 
-0.62 0 0 
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5.11 0 0 

» S 2 10 8 0 0.16 0 

0 0 0.39 

The variance of the angLe x is 

cr 2 = 2.26 X 10 4 s . 
x 

Withthe current accuracy of s = 3.0 x 10“ 6 [10] 

or 2 = 6,78 x 10 -2 rad = 118km/1738km 
x 

i.e. , a standard deviation of about 120 km of arc on the lunar surface. 
tv. r>Kf-air. on or.o. ira cy of 1km, potential coefficients would need to be 


2.26 x 10‘ 


lkm 10” 4 
1738 km 2.26 


= 2.5 X 10" 


i.e. , an improvement of the current aecuracyby a factor of 100 would be 
needed. 
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